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ABSTRACT 

The GJ 436 system contains a transiting planet GJ 436 b which is a hot analogue of Neptune on 
an eccentric orbit. Recently, two additional transiting sub-Earth planets have been postulated in the 
literature. We observed three transits of GJ 436 b over the course of 3 years using two-meter class 
telescopes, each with a photometric precision better than one millimagnitude. We studied system 
dynamics based on the existence of the additional planets. We redetermined system parameters, 
which were in agreement with those found in the literature. We refined the orbital period of GJ 436 b 
and found no evidence of transit timing variations. The orbital motion of the GJ 436 c planet candidate 
was found to be significantly affected by the planet b with variations in transit times at a level of 20 
minutes. As the orbital period of the GJ 436 d planet candidate remains unknown, our numerical 
experiments rule out orbits in low-order resonances with GJ 436 b. The GJ 436 system with the hot 
Neptune and additional two Earth-like planets, if confirmed, would be an important laboratory for 
studies of formation and evolution of planetary systems. 

Key words: planetary systems - stars: individual: GJ 436 — planets and satellites: individual: 
GJ 436 b. GJ 436 c, GJ 436 d 


1. Introduction 

The GJ 436 planetary system, with one confirmed planet and two proposed 
planets, may be the first multiple-transiting-planet system initially discovered from 
the ground. The host star - an M3 dwarf with an age of several Gyr, located 10 
pc from the Sun - was found to have a Neptune-mass planetary companion on a 
2.64 day orbit with the precise radial velocity (RV) technique (Butler et al. 2004). 
Additional Doppler measurements found an orbit that was far from circular, with 
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an eccentricity of Cb = 0.16 ±0.02 (Maness et al. 2007) and the semi-major axis of 
which is Qb = 0.0285 AU. Gillon et al. (2007b) later detected planetary transits in 
photometry of the host star. Combining the spectroscopic and photometric data al¬ 
lowed those authors to determine the planet’s mass of Mb = 22.6 ±1.9 (Earth 
masses) and radius of Rb = 25200 ±2200 km = 3.95 ±0.34 (Earth radii). 
These properties showed that GJ 436 b must be a planet similar to Uranus or Nep¬ 
tune, but much closer to its host star. Due to this proximity, the planet’s atmosphere 
is hot, with an equilibrium temperature between 520 and 620 K. These features of¬ 
fer unique opportunities for a number of follow-up observations conducted with 
ground-based and space-born facilities. 

Maness et al. (2007) and Deming et al. (2007) have postulated that the non¬ 
zero eccentricity could be caused by an unseen planetary companion. Based on 
available RV data. Ribas et al. (2008a) found a sign of additional planet close to 
the 2:1 mean motion resonance with GJ 436 b. Coughlin et al. (2008) indicated 
that orbital inclination, transit depth, and transit duration may exhibit variations 
excited by gravitational influence of another planet in a non-resonant orbit. On 
the other hand, dynamical studies of Alonso et al. (2008) and Bean & Seifahrt 
(2008) eliminated the proposed 2:1 configuration, and placed physical limitations 
on possible configurations of the second planet in the system. 

Using the Spitzer Space Telescope, Stevenson et al. (2012) detected addtional 
shallow transit-like features in the light curve of GJ 436. They propose that these 
features could be caused by transits of two additional planets with radii of 0.7 
/?0. The orbital period of the planet candidate GJ 436 c (originally labelled UCE- 
1.01) was found to be 1.37 d, while the orbital period of GJ 436 d (UCE-1.02) could 
not be determined because only two transits were observed. Assuming a range of 
bulk densities typical for terrestrial planets (/.e,between 3 and 8 g cm^^), the 
masses of both planet candidates were constrained to 0.15-0.40 - much below 

a detection threshold of current most advanced RV surveys. The orbital period 
ratio of GJ 436 b and GJ 436 c is 1.94, close to a 2:1 orbital resonance which is 
generally thought to produce strong transit time variations (TTVs). Recent studies 
by Eanotte et al. (2014) and Stevenson et al. (2014) show that the existence of both 
planet candidates still remains disputable. 

It is believed that future space-based instruments will provide opportunities to 
confirm the multi-planetary architecture of the GJ 436 system. We note, however, 
that TTV observations from the ground and dynamical studies may place some 
constraints on possible planetary configurations. In this paper, we explore the pos¬ 
sibility of the proposed multi-planetary architecture, using new upper limits on the 
TTVs of GJ 436 b based on our photometric observations and combining this with 
a dynamical model of the system. 
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2, Observations and data reduction 


We observed the 2011 Jan 04 transit of GJ 436 b with the 2.2-m teleseope at 
the Calar Alto Observatory (Spain) as a baek-up target within the program FI 1-2.2- 
008, whieh was foeused on transit timing of the WASP-12 b planet. The Calar Alto 
Faint Objeet Speetrograph (CAFOS) in imaging mode was used to aequire the light 
eurve in the Cousin R filter. To shorten the read-out time, the original field of view 
(FoV) was windowed fo 5.3' x 6.5'. Observafions were oeeasionally affeeled by 
passing fhin elouds. We observed a seeond fransif on 2012 Feb 02 wifh fhe Nordie 
Opfieal Teleseope (NOT) af fhe Observaforio del Roque de los Muehaehos, La 
Palma (Spain), as a baekup large! of fhe P44-102 observing program (OPTICON 
2011B/003), fhe goal of whieh was aequiring high-preeision fransif lighl eurves 
for WASP-12 b’s Iransils. The Andalueia Fainl Objeel Speelrograph and Cam¬ 
era (ALFOSC) was used in an imaging mode, windowed fo fhe 6.4' x 6.4' FoV. 
Phofomelrie lime series was laken in fhe Bessel R filler under exeellenl wealher 
eondilions. 

Bolh inslrumenls used 2x2 binning for faster readoul. Auloguiding guar- 
anfeed fhal stellar eenlroids did nol ehange Iheir loealions on fhe deleelor malrix 
during eaeh run. We defoeused fhe leleseopes signifieanlly, erealing doughnul- 
like stellar profiles and spreading fhe lighl over many pixels, whieh enabled longer 
exposure limes and redueed fhe impael of llal-fielding imperfeelions. This also im¬ 
proves duly eyele, as fhe ratio of time during exposure versus read-oul is improved, 
allowing more photons to be galhered during fransif and Ihereby inereasing fhe pho¬ 
fomelrie effieieney. We redueed fhe observations wifh slandard proeedures, inelud¬ 
ing de-biasing and fial-fielding using sky fials. We performed differenlial aperlure 
pholomelry wifh respeel fo nearby eomparison slars BD-i-27 2046 and TYC 1984- 
1884-1. Radii of fhe aperlure and baekground ring were empirieally optimized to 
produee fhe smallesl sealler in fhe oul-of-lransil lighl eurve. 

We oblained an additional fransif lighl eurve from observafions made on 2014 
Mar 23 from low-resolution speelra aequired wifh fhe NOT/ALFOSC. For Ihese 
observations we ehose fhe 2x2 binning mode, a readoul speed of 200 pixeFsee 
wifh a gain of 0.327 e^/ADU and a readoul noise of 4.2 e^/pixel. We used AL¬ 
FOSC in ils long-slil speelroseopie mode, seleeling fhe grism #4 whieh eovers fhe 
speelral range 3200 - 9100 A and a euslom-buill slil of a widlh of 40 aresee. We 
ehose BD-i-27 2046 as referenee, whieh has similar brighlness and is loealed al a 
disfanee of 3.81 aremin from GJ 436. The position angle of fhe referenee slar wifh 
respeel fo fhe largel was equal to 37.39°. Observafions began al 01:31 UT (35 min- 
ules before ingress) and ended al 03:18 UT (10 minules afler egress). The exposure 
time was sel to 60 seeonds and fhe readoul time of inslrumenl was 9 s, meaning 
we eolleefed approximately one speelrum every 69 s. The dafa reduelion (bias and 
flal-field eorreelions, exlraelion of fhe speelra and eorresponding ealibralion ares 
as well as wavelenglh ealibralion using Ihe He and Ne lamps) was made using an 
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IRAF0 script written for NOT/ALFOSC long-slit data. Optimal extraction of the 
spectra was obtained using an aperture width of ± 10 binned pixels, which corre¬ 
sponds to 7.6 arcsec on the detector. This is 2 to 5.5 times the raw seeing during 
the observations (1.4 - 3.7 arcsec). The light curve was constructed using 1/3 of 
spectra of the target and reference star centered at the maximum of the GJ 436 
spectrum. The maximum of the spectral energy distribution was found at 790 nm 
that corresponds to the central passband of the photometric I filter. 

Differential atmospheric extinction and differences in spectral types of the tar¬ 
get and comparison stars, as well as instrumental effects caused by the field dero- 
fafor in NOT dafa, are expecfed fo produce phofomefric frends, whose time scale 
is similar fo fhe durafion of fhe run. The de-frending procedure was done wifh fhe 
JKTEBOP code (Soufhworfh et al. 2004a,b) by tiffing a second-order polynomial 
funclion of lime along wifh a frial fransif model and fhen subfracfing fhe resulfing 
polynomial from fhe lighf curve. Magnifudes were Iransformed info fluxes, which 
were normalized fo unify oulside of fhe fransif. The fimesfamps in geocenlric Ju¬ 
lian dales in coordinafed universal time (UTC) were provided by a GPS syslem 
and verified wifh fhe nelwork lime prolocol soflware, fhen converted fo barycenfric 
Julian dales in barycenfric dynamical time (BJDTDB^ Easlman et al. 2010). The 
observing delails and dafa qualify characlerislics are given in Table 1. 


Table 1 

New transit light curves reported for GJ 436 b: date UT is given for the middle of the transit, X 
shows a course of changes in airmass during transit observations, texp is the exposure time, F is 
the median number of exposures per minute, pnr is the photometric scatter in millimagnitudes per 
minute of observation, as defined by Fulton el al. (2011). 


# 

Date UT 

Telescope 

X 

texp (s) 

r 

pnr (mmag) 

1 

2011 Jan 04 

2.2-m Calar Alto 

1.08^ 1.02-^ 1.07 

12 

1.76 

0.6 

2 

2012 Feb 02 

2.6-m NOT 

1.08 1.001.06 

20 

2.07 

0.3 

3 

2014 Mar 23 

2.6-m NOT 

1.01 ^ 1.18 

60 

0.87 

0.7 


3. Results 

We used different codes to model the transit light curves and analyze the RV and 
timing datasets. We obtained the final resulls of Ihis sludy by iferafing fhe differenl 
codes unfil we reached convergence. Sections 3.1 and 3.2 confain a description of 
fhe sleps involved in each iferafion. 

*IRAF is distributed by the National Optical Astronomy Observatory, which is operated by the 
Association of Universities for Research in Astronomy (AURA) under cooperative agreement with 
the National Science Foundation. 






4 


3.1. Transit light curves 

We modeled the new sub-millimagnitude precision light curves simultaneously 
with the Transit Analysis Packag^ (TAP, Gazak et al. 2012). This code em¬ 
ploys the Markov Chain Monte Carlo (MCMC) method, including the Metropolis- 
Hastings algorithm and a Gibbs sampler, to find the best-fit parameters of a transit 
light curve approximated by the analytical model of Mandel & Agol (2002). The 
time-correlated noise in data (so called red noise) is investigated with the Carter & 
Winn (2009) wavelet parametrization. This approach is known to provide conserva¬ 
tive uncertainty estimates. TAP parametrizes the flux distribution across the stellar 
disk with a quadratic limb darkening (LD) law (Kopal 1950). The values of linear 
and quadratic LD coefficients, u\ and U 2 respectively, were linearly interpolated 
from tables of Claret & Bloemen (2011) with an on-line tooJl of the EXOFAST 
package (Eastman et al. 2013). The stellar parameters for GJ 436 were taken from 
von Braun et al. (2012), assuming solar metallicity. 

The fitting procedure kept the orbital inclination ij ,, the semimajor-axis scaled 
by stellar radius ah/R*, and the planetary to stellar radii ratio Rh //?* as free param¬ 
eters, linked together for all light curves. The mid-transit times were determined 
independently for each light curve. The orbital period was fixed af a value from a 
refined ephemeris. The ED coefficienfs were allowed fo vary around fhe fheorefical 
values under fhe Gaussian penally of a = 0.05, independenlly for all fhree dafasefs. 
This allowed us fo accounf nof only for differences belween R and I bands, buf also 
for any possible differences in insfrumenfal implemenfafion of /?-band fillers. The 
orbifal eccenfricily et and longifude of periaslron (at, were faken from fhe dynam¬ 
ical model discussed in Seel. 3.2. The new lighl curves wilh fhe besl-fil Iransil 
model are plolled in Fig. 1. 

Our new Iransil observations prolong fhe fimespan covered by observations and 
allow fhe fransif ephemeris fo be refined. The derived mid-lransif limes were com¬ 
bined wilh fhe lileralure ones fo calculale new reference epoch Tq and orbifal period 
Pb. 

The paramefers of our besl-fif model are given in Table 2, logefher wilh recenf 
lileralure resulls for comparison. We also lisl fhe resulfing system properties: fhe 
Iransil parameter bb , defined as 


bh 


at 1 - 

- - -cos lb , 

7?* l+eb cos (Ob 


( 1 ) 


and fhe mean stellar density p*, which can be direcfly calculaled from Iransil ob¬ 
servable properlies wilh a formula derived from Kepler’s fhird law 


p* — 


371 



( 2 ) 


^http://ifa.hawaii.edu/users/zgazak/IfA/TAP.html 
^ http:// astroutils. astronomy.ohio- state, edu/ exofast/limbdark. shtml 
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Fig. 1. New transit light curves for GJ 436 b with the best-fit model, plotted with solid lines. The 
residuals are shown in bottom parts of the panels. 

where G is the gravitational eonstant. 


Table2 

Parameters of the GJ 436 system derived from modeling transit light curves. 


Parameter 

This work 

von Braun et al. (2012) 

Lanotte et al. (2014) 

ib n 

Ob/R* 

16 

1 -T '70+0.46 
1D./D_o43 

86.6+1 

86:858™ 

14.54+T 

Rb/R* 

bb 

P* (po) 

U. /JO_o,041 

^■^'-0.47 

0.0833 ±0.0002 

n oco+O-OOS 

U.QJD_qqq 3 
c q'7+0.30 
^•■^'-0.27 

0 7079 + 0*3053 

C Q1 +0.17 
^•^^-0.18 

To (BJDxdb) 

2454510.80162±0.00007 

2454510.80096±0.00005 

- 

Pb(d) 

2.64389754 ±0.00000043 

2.t)4389520_Q_QQQQQQ5g 

9 f^J.98Q8n'^+0 00000027 


3.2. Dynamical model with the GJ 436 c planet candidate 

The linear fit to mid-transit times for GJ 436 b results in a reduced equal 
to 5.5, which indicates a marked departure from a linear ephemeris (Stevenson et 
al. 2012). The periodogram analysis reveals no statistically significant periodic 
signal, so stellar activity, systematic effects, or underestimated timing errors may 
be a source of a spurious timing variations. An upper limit on the amplitude of any 
periodic signal in transit timing was found to be 0.0002 d. 

Mid-transit times for the GJ 436 c planet candidate show noticeable departure 
from a linear ephemeris (Stevenson et al. 2012) that suggests planet’s orbital motion 
is perturbed by the planet b. To study the mutual interactions between both planets, 
a two-planet dynamical model was constructed with the Systemic code in version 
2.16 (Meschiari et al. 2009). We used 113 RV measurements from Knutson et al. 
(2014), acquired with the HIRES echelle spectrometer at the Keck I telescope be¬ 
tween January 2000 and December 2010. We also used 171 HARPS RVs obtained 
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between January 2009 and April 2010 from Lanotte et al. (2014). For GJ 436 b, 
we used mid-transit times for 37 epoehs available in the literature (Gillon et al. 
2007a,b; Shporer et al. 2009; Caeeres et al. 2009; Deming et al. 2007; Alonso 
et al. 2008; Knutson et al. 2011; Pont et al. 2009; Bean et al. 2008; Ribas et 
al. 2008b; Coughlin et al. 2008; Ballard et al. 2010; Beaulieu et al. 2011) and 3 
new determinations reported in this paper. In addition, we used 15 mid-oeeultation 
times from Stevenson et al. (2010) and Knutson et al. (2011). The best-fit Newto¬ 
nian solution was found with a differential evolution algorithm using 5000 steps, 
followed by a number of iterations of Levenberg-Marquardt optimization. The 
Runge-Kutta-Fehlberg (RK45) algorithm was used to integrate equations of mo¬ 
tion with an aeeuraey requirement of 10^^^. The bootstrap method with 10^ trials 
was used to estimate parameter uneertainties, ealeulated as median absolute de¬ 
viations. The orbital periods, eeeentrieities, and arguments of periapsis were left 
free for both planets. The mass of the planet eandidate was fixed af fhe value of 
0.28 (Earfh masses) faken from Stevenson et al. (2012). The mass of plane! 
b was allowed fo vary. Orbifal inelinafions were faken from our fransif lighf-eurve 
analysis (Seef. 3.1) for GJ 436 b and from Stevenson et al. (2012) for fhe planel 
eandidate. The paramefers of fhe besl-fil dynamieal model are listed in Table 3. 


Table3 

Orbital parameters for the GJ 436 b planet and GJ 436 c planet candidate from the two-planet 
dynamical model. The values are given for the epoch JD 2455959. 


Parameter 

GJ 436 b 

GJ 436 c 

Orbital period (d) 

Semi-major axis (AU) 

Orbital eccentricity 

Longitude of periastron (deg) 
RV amplitude (m/s) 

Mass (M 0 ) 

2.64388312 ± 0.00000057 

0.0291 ±0.0015 

0.13827 ±0.00018 

351.00 ±0.03 

17.09 ±0.22 

22.1 ±2.3 

1.365960 ±0.000012 

n 01S71 +0.00097 
O-Olo/1-0.00094 

0.1166 ±0.0046 

36.50±0.41 

0.3* 

0.28** 

* predicted value 

** value taken from Stevenson et al. (2012). 


The SWIFT’s Regularized Mixed Variable Sympleefie (RMVS) infegrafor was 
used fo fraee fhe evolufion of orbifal paramefers, whieh exhibif oseillafory pafferns 
as a funefion of time. The besf-fif dynamical model is sfable in a fimescale of 10^ 
yr (over 130 million orbifs of GJ 436 b) and yields sfrong consfrainfs on fhe ec- 
cenfricify of fhe GJ 436 c planel candidale. Ifs besl-fil value is similar fo lhal one 
for GJ 436 b bul if is expecled fo oscillate belween a value as small as 0.02 and 
0.19 wilh a period of 35.4 yr (Fig. 2a), af anti-phase wilh marginal variations in eb 
(belween 0.137 and 0.139; fhe range is smaller by a factor of fhe mass ratio lhaf is 
a consequence of conservation of momenlum). The predicted variation in orbifal 
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2020 2040 2060 2080 2100 

Time (yr) 

Fig. 2. Evolution of the orbital eccentricity and inclination for the GJ 436 c planet candidate (panels 
a and b) and difference in arguments of periapsis between GJ 436 b and c. 

inclination for GJ 436 c is 2.6° (Fig. 2b). The inclination of GJ 436 b’s orbit was 
found to decrease with a rate of 0.025 degree per century - far under the detec¬ 
tion threshold of current transit observations. The ratio between the orbital periods 
of both planets, Pb/Pc = 1-94, suggests the planets are close to a 2:1 commensii- 
rability and could be trapped in a mean motion resonance. The evolution of the 
difference between arguments of periastron, defined as Aco = 00^ — COc GS plotted in 
Fig. 2c. The periastrons were found to be in apsidal alignment around 0° and li- 
brate with an amplitude of 56°. Such conditions are generated by the linear secular 
coupling and prevent both planets from close encounters which could destabilize 
the system. The eccentricity-type resonant angles, defined as a linear combinafion 
of mean longifudes X and argumenfs of periasfron, 

^b = Xc — 2Xb + (Ob (3) 

9c = “ 2Xb -F (Oc (4) 

show no librafion, so bofh planefs are nol in a dynamical resonance. 

Figure 3 shows fransif liming residuals for GJ 436 b and fhe GJ 436 c planel 
candidale, produced by mufual gravifafional inleraclions. For bofh planefs, fhe 
signal is periodic wifh a period of 41.2 d. The range of variation is subsfanlial for 
fhe GJ 436 c planel candidale wifh a value of almost 20 min. Timing of GJ 436 b is 
predicted to be modulated with an amplitude of 16 s - a signal which is under the 





5200 5400 5600 5800 

BJDtdb-2450000 


Fig. 3. Timing residuals for transits of the GJ 436 b planet and the GJ 436 c planet candidate, resulting 
from mutual gravitational interactions. 

detection threshold of the current timing dataset. 

3.3. Constraints on the orbital period of the GJ 436 d planet candidate 

The period of time between the two observations of the transits of candidate 
GJ 436 d, Xd = 151.5 d, must be a multiple of a true orbital period Pd. Stevenson 
et al. (2012) estimate the upper limit for Pd to be 5.56 d, based on 1-a uncertainty 
in transit duration of the planet candidate. We note however, that the fact that 
the orbital inclination is unknown weakens their argument. Periods shorter than 
Pb may be excluded because no transit signature was found in the 2.4-day long 
continuous light curve of GJ 436 acquired with the Spitzer Space Telescope at 
8 qm (Stevenson et al. 2012). We used GJ 436 b’s transit timing and system’s 
dynamical stability to put constraints on Pd. 

Using the Systemic code, we conducted a numerical experiment in which the 
planet d was inserted into the two-planet system derived in Sect. 3.2. The mass of 
the planet d was set to 0.27 as given by Stevenson et al. (2012). The orbital 
eccentricity, ed , varied between 0.0 and 0.3 with a step of 0.05. The initial value 
of the argument of periastron was set equal to co^ and initial orbital longitude was 
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shifted by 180° with respect to the value for GJ 436 b at the epoch 0. The orbital 
inclination of the planet d was the same as for GJ 436 b. The value of Pd satisfied 
the relation 

Xd = n-Pd, (5) 

where n is 1, 2, 3... Each configuration was integrated with the RK45 algorithm, 
covering 2500 days, /.e,the time span of the transit observations for GJ 436 b. 
The times of mid-transit times were extracted and the amplitude of periodic devia¬ 
tions from a linear ephemeris was calculated. In addition, each configuration was 
checked for the dynamical stability in 10^ yr with the SWIFT’s RMVS integrator. 
Configurations which were found to be unstable, mainly with Pd < 3/2Pb, were 
skipped. 

The exemplary results for Cd equal to 0.05, 0.15, and 0.25 are shown in Fig. 4. 
Configurations with Pd/Pb close to 3:2 and 2:1 resonances would generate TTV 
signals above a detection threshold even for low-eccentricity orbits of the planet d. 
For greater eccentricities, configurations with resonances close to 3:1 or 5:2 would 
produce detectable TTVs. 


4. Concluding discussion 

Our new observations of GJ 436 b’s transits allowed us to refine transit ephemeris 
and to redetermine system parameters. They prolong the timespan covered by ob¬ 
servations by a factor of 3. They were combined with photometric and Doppler 
data from the literature in order to study dynamics of the system with two low- 
mass planet candidates. 

The interferometric measurements for the GJ 436 host star give a precise value 
of stellar radius of = 0.455 ±0.018 Rq (von Braun et al. 2012), and when com¬ 
bined with a mean stellar density of p* = 4.971 q 47 p© results in a stellar mass of 
M* = 0.47 ±0.07 Mq . This value is between 0.452 +qq}2 Mq and 0.507+QQg2 
reported by Torres (2007) and von Braun et al. (2012), respectively. It is also still 
consistent within a 1-a range with the value of 0.556('(QQg5 M© given by Fan- 
otte et al. (2014). The deduced stellar surface gravity is logg* = 4.792 ^qq 44 in 
cgs units, consistent with logg* = 4.83 ±0.03 derived by von Braun et al. (2012). 
Our determinations of the GJ 436 b’s radius Rb = 0.372 ±0.015 Rjup and mean 
density pb = 1-35 ±0.22 pjup agree with the recent literature values within 1-a 
(0.369 ±0.015 Rjup and 1.551 o}q reported by von Braun et al. 2012, and 
0.366 ±0.014 Rjup and 1.6 pjup given by Fanotte et a/. 2014). Our value of plan¬ 
etary mass Mb = 22.1 ±2.3 M© is in a perfect agreement with Mb = 22.6 ±1.9 
M© given by Gillon et al. (2007b), but seems to be slightly underestimated com¬ 
paring to 24.8^2 5 of Braun et al. (2012) and 25.4^2 0 of Fanotte et 
al. (2014). This is a direct consequence of higher stellar mass determined in both 
studies. 
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10 100 


Orbital period (d) 


Fig. 4. Amplitude of a hypothetical TTV signal for GJ 436 b induced by the GJ 436 d planet candidate 
as a function of possible values of for selected values of . Only configurations that are stable 
for a time scale of 10^ yr are shown. The horizontal dashed lines show a detection threshold based 
on observations. Vertical dotted lines mark selected to Pi, resonances. 


Mid-transit times reported in this study follow a linear ephemeris with residuals 
smaller than 2a. The transit from 2011 Jan 04 (JD 2455565.7), the most outlying 
from the linear ephemeris in our sample, eould be affeeted by stellar aetivity at the 
ingress phase beeause photometrie residuals from the transit model seem to exhibit 
some distortion. The laek of a periodie or semi-periodic TTV signal indicates that 
GJ 436 b is not noticeably perturbed by gravitational interactions with other bodies 
in the system. On the other hand, the orbital motion of the GJ 436 c planet can¬ 
didate is significantly influenced by GJ 436 b. Our numerical model predicts that 
periastrons of both planets are in apsidal alignment. The TTV signal for GJ 436 c 
is expected to be at the level of 20 min and constrains the orbital eccentricity of 
the planet. Further precise transit observations for GJ 436 c will shed new light for 
system dynamics. 

The orbital period of the GJ 436 d planet candidate remains unknown but the 
lack of the TTV signal for GJ 436 b put some constraints on it. Configurations 
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in which both bodies are close to low-order resonances are unlikely because they 
would generate detectable TTVs for GJ 436 b. A system with planets out of reso¬ 
nance is in line with statistical studies based on data from the space-based Kepler 
transit survey, which show that configurations with planets in close proximity to a 
resonance are not favoured (Fabrycky et al. 2014). 
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